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1. Introductory Comments 

This document reports the work completed on the NASA research grant NAG-1-612, 
entitled "Response of Joint-Dominated Space Structures." The work falls into two 
categories: 1) developing an efficient method for calculating the transient response 

of a nonlinear system such as a large joint-dominated space structure and 2) 
investigating the effect of gravitational loading and joint scaling on the dynamic 
response of model structures. The results of the investigations on these two topics 
are reported in Section 2 and Section 3, respectively. 


2. Transient Response of Joint Dominated Space 
Structures: A New Linearization Technique 


G. A. Foelsche 
J. H. Griffin 
J. Bielak 


Carnegie Institute of Technology 
Carnegie Mellon University 
Pittsburgh. PA 15213 


Transient Response of Joint Dominated Space 
Structures: A New Linearization Technique 


by 

G. A. Foelsche, J. H. Griffin, and J. Bielak 
Carnegie Mellon University 
Pittsburgh, Pennsylvania 15213 


Abstract 

A new, efficient linearization method is presented for calculating the transient 
response of nonlinear systems due to initial disturbances. The method is an 
extension of the describing function approach in which the steady-state response of 
the system is calculated by representing the nonlinear element, typically joints in the 
the case of space structures, by impedances which are functions of the amplitude of 
response. Thus, the problem of solving the differential equation for the steady-state 
response becomes one of solving a set of nonlinear algebraic equations involving the 
steady-state amplitudes and phases of the system. It is shown that for the transient 
case the steady-state impedances can be averaged over the range of response in 
order to provide equivalent values of stiffness and damping that, for a given set of 
initial displacements, may be treated as being constant for purposes of calculating 
system response. 

Single-degree-of-freedom systems are used first to demonstrate the method and 
then to develop an approach for optimizing the joint's characteristics so as to 
minimize transient response times. The use of this method for response estimation 
and optimization in multiple-degree-of-freedom systems is investigated subsequently. 
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stiffness 

stiffness of nonlinear friction element 

effective value of instantaneous system stiffness 

mass 

normal load on friction contact 
time 

half period 
eigenvectors 
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fraction of critical damping 

nondimensionalized quantities 

quantities averaged over a range of response 
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Introduction 

An important problem related to the design of space structures is that of predicting 
their dynamic response. This process is particularly difficult when dealing with 

structures that contain a large number of joints that exhibit nonlinear, hysteretic 
behavior. For example, this may be the case for prefabricated truss structures that 
are designed to collapse into a dense package for transportation to orbit. The truss 
is then expanded in space by utilizing joints that are especially designed to rotate 
and lock into place. The dynamic response of such structures is said to be joint 
dominated if the amount of damping or the stiffness of the system is strongly 
affected by the joints' behavior. If damping in the system is primarily due to joint 
hysteresis then joint behavior controls the amplitude of the steady state response as 
well as the rate at which transients decay. Additionally, in some cases, joint 
flexibility can significantly reduce the stiffness of the structure, thus reducing its 
natural frequencies and altering the associated mode shapes. This paper discusses a 
new approach that may be used to efficiently estimate the transient response of such 
systems. 

The transient response of nonlinear systems is usually calculated by time 
integration methods that employ finite differences in time; see, for example, Hughes 
[1]. This approach has two disadvantages when it is applied to the design of joint 
dominated structures. The first is that it is computationally intensive. This would 
be especially the case for the type of complex three-dimensional truss structures that 
are proposed for space applications since they have a large number degrees of 
freedom. Secondly, the problem is nonlinear and, consequently, the solutions lack 
generality. For example, the rate of decay of a transient would depend on the 
specific magnitude and distribution of the assumed initial displacements and 
velocities of the structure. Since the number of degrees of freedom is large, it is 
not reasonable to consider all possible initial conditions. So one difficulty faced by 
the design engineer is in trying to select those conditions that are of critical 
importance to the design and to simulate the corresponding system response when 
there is such a large number of possible cases to consider. Neither the selection nor 
the simulation process is particularly feasible if time integration Is the only 
procedure available for calculating the response of the system. 

In this paper an approximate method is developed for estimating the transient 
response of nonlinear systems in terms of linearized modes of response. Its 

advantages are that it is computationally more efficient than the time integration 
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method and that it is possible to view the design problem in the more traditional 
physical terms of modal response. For example, if it is most important to damp the 
larger amplitude, low frequency response one can easily focus on that issue by 
isolating the response of the first few linearized modes. Consequently, by utilizing 
this approach the design problem should become more tractable. The major drawback 
of the approximate method is loss of accuracy. It is our view that both approximate 
methods and time integration have their roles in design. Approximate methods 
provide efficient tools for performing parametric studies and they supply physical 
insights into how to optimize system performance that are not easily inferred from 
strictly numerical methods. Time integration provides a method for assessing the 
accuracy of the approximate solution for key simulations and for fine tuning the final 
design. 

In the procedure presented here the nonlinear system is approximated by an 
equivalent linear system in which the system parameters are constant over the^ range 
of transient response. The method is an extension of the describing function 
approach used to calculate the steady-state harmonic response of nonlinear systems. 
In the describing function approach the response is assumed to be essentially 

harmonic and the nonlinear element is represented by impedances which depend on 
the amplitude of response. As a result, the problem of solving the differential 
equation for the steady-state response becomes one of solving a set of nonlinear 
algebraic equations for the steady-state amplitudes and phases associated with the 
various degrees of freedom of the system. In the transient case considered here the 

steady-state impedances are averaged over the range of response in order to provide 

equivalent values of stiffness and damping that, for a given set of initial conditions, 
may be treated as constants for purposes of calculating system response. We refer 
to this approach as the Amplitude Averaging (AA) Method. Once equivalent 

parameters are identified for the system, conventional methods can be employed for 
analyzing the resulting linear system. Related studies have been summarized, for 
instance, by Iwan and Gates [2]. 

The AA Method is derived from an efficient time integration procedure presented 
by Sinha and Griffin [3] in which single time steps were used to step from one 
peak of the oscillation to the next. Their approach in turn, was based on the idea 
that the response may be approximated as a sinusoid in which the amplitude and 
phase vary slowly with time (see, for example, Caughey [4]). 


In the first part of this paper the AA Method is illustrated by applying it to a 
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single-degree-of-freedom (SDF) system exhibiting the bilinear hysteretic behavior 
typically associated with Coulomb friction. While the method is not restricted to 
this type of nonlinearity this behavior was selected for analysis because it is 
representative of the type of severe nonlinearity that occurs in actual joints. A half 
cycle method similar to that used by Sinha and Griffin is presented in order to 
illustrate the linearization process and in order to develop "instantaneous" values of 
the nonlinear element's stiffness and damping. The instantaneous values of the 
element's parameters, which effectively characterize the joint's properties at a given 
amplitude of response, are then averaged over a range of amplitudes in order to 
calculate the constant stiffness and damping values used in the AA Method . 1 Since 
an equivalent, constant damping has been determined for the system it may be used 
to select joint characteristics (the friction slip load in this example) so as to 
maximize average joint damping and minimize transient response times. In the 

second part of the paper a general approach that may be used for multiple-degree-of- 
freedom (MDF) systems is given and applied to the two body problem. In each case, 
the accuracy of the approach is assessed by comparing results from the approximate 
method with those obtained using standard time integration methods. 

Single Degree of Freedom Systems (SDF) 

In order to demonstrate the AA Method, we develop a solution for the bilinear 
hysteretic SDF system depicted in Figure 1. The equation of motion for the system 
is 


mx ♦ cx ♦ kx « - f (1) 

ft 

where f is the nonlinear force from the friction element and the dots represent 
differentiation with respect to time. During oscillation, the friction joint remains 
locked until the magnitude of the spring force, |k d (x-y)|, equals the friction force 
fiN. The joint then slips with a constant resistive force of magnitude /iN until the 
mass reaches an extremum of oscillation, at which point the joint locks up again. 
The magnitude of the relative displacement (x-y) required to cause slip is designated 
as x , where 


\n the example used here to demonstrate the *>proach, the instantaneous joint properties are calculated 
analytically. Actual joint properties may be calculated from laboratory tests that measure the joint's steady state 
hysteresis curves. See Crawley [5] for an example of steady state joint characterization. 
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The Half-Cycle Method 

If the transient response is approximately a sinusoid that has an amplitude and 
phase which vary slowly with time, then over a limited time span it may be 
approximated as x « B cos0, where 6 * »t - f. We assume that the nonlinear force 
f exhibits the same periodicity. Expressing f in a Fourier series, we obtain (see 

n n 

Menq and Griffin, [6]) 

f « F (B) cos# * F (B) sin0 + (higher harmonics) (3) 

n c s 

where, for the case of Coulomb friction considered; 


F (B) * f 0*(B) - ]■ sin20*(B) 1 
c . L 2 J 


(4) 


4 M N . f, N v 

F (B) * — — ( 1 - ) 

* n ' k£ ' 

0 


(5) 


0*{B) • cos” 1 ( 1 - ) 

\ k d B / 


0 £ 6 S, v 


( 6 ) 


In the case of experimental joint data, the Fourier coefficients can be calculated 
numerically from hysteresis curves using numerical integration. 

If we keep only the fundamental harmonics and truncate the expression (3) after the 
first terms, we can express f as 


f 

n 



-F (B) . 
x + — - — x 
u B 


(7) 


Thus (1) becomes. 
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mx + 



x a 0 


( 8 ) 


or. 

mx + c eff (B)x + k eff (B)x = 0 (9) 

where c .. and k .. are defined in a manner consistent with (8) and because they 
eff eff ' 

pertain to a specific amplitude are referred to as the instantaneous damping and 
stiffness of the bilinear spring. Alternatively, we may write (9) as 

x + 2{w x + w 2 x = 0 (10) 

n n 

where £ and w are both functions of B 

n 


£(B) 


2(m k <B)) W 
eff 


(id 


(B). 

n ' m • 


( 12 ) 


The damped natural frequency of this system may be similarly defined as 


» d (B) * 


M - t 2 (B))k j „(BI . „ 2 
' m • 


(13) 


It is observed that during transient oscillation of the nonlinear system, the response 
is similar to the decaying sinusoid seen in linear analysis. It is reasonable to 
assume that the motion of the nonlinear system from one extremum of oscillation to 
the next extremum is representable as the decay of a linear system over a half-cycle. 
As an example, we consider a system decaying from initial conditions of some 
initial displacement B q and zero initial velocity. Let B. denote the amplitude of the 
i^ extremum (occuring at time t.) and B. + 1 denote the amplitude of the next extremum 
(at time t j+1 ). Then from linear theory, the time elapsed between one extremum and 
the next is approximately 
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At = t 


i+1 


- t. 
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n 

w (B.) 

d i 


(14) 


and the relationship between successive peaks is 

-c »/<i - c.V ' 2 

B ,., - - B , • ' ' 1151 

Given the initial amplitude B , we can efficiently estimate successive extrema and 
the time increments at which they occur from (14) and (15). Since the extrema occur 
each half cycle, we refer to this method as the Half-Cycle Method. 


This Half-Cycle Method is an efficient way of approximating a numerical time 
integration to find extrema of transient oscillation. Numerical simulations of SDF 
systems show that this method is accurate for nonlinear systems in which t = k rf /(k + 
kj £ 0.5 and provides a reasonable approximation for the amplitude of response for 
t > 0.5. The Half-Cycle Method results were compared to more accurate solutions 
generated by fourth order Runge-Kutta time integration. The comparison of the half 
cycle estimates to the numerically generated "exact" solutions is shown in Figures 
2a-d where the Half-Cycle estimates of the extrema of response are shown as points 
(the exponential decay envelopes pictured are a result of the new Amplitude 
Averaging Method and will be discussed later). For nonlinearities of < > 0.5, the 
system experiences an offset that is not accounted for in the Half-Cycle Method. 
However, the peak-to-peak amplitude estimates are approximated reasonably well 
despite this offset. It may be observed that for these simulations this approach 
yields conservative estimates of system behavior in that it overestimates the 
amplitudes of response. 


Amplitude Averaging Method 

We now introduce a new linearization called the Amplitude Averaging (AA) Method. 
In this approach, the half-cycle values of damping and stiffness are averaged over 
the entire response range of interest. Consequently, while these average values are 
nonlinear functions of the initial displacement, they are constants as far as the 
transient response analysis is concerned. 

The AA Method is also based on equations (3) through (6). Linear, constant 
parameters are derived by averaging c #ff and k 0ff from (9) over the nonlinear range 
of oscillation, i.e.. 
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eff o e o 


f iff 

J Bar 

crit 


dB/(B - x ) 

o crit 
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(B) 
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g o crit 


crit 


(16) 


(17) 


where w in (16) is given by (12). After averaging, (8) becomes 

mx + (c + c (B ))x + (k + k (B ))x * 0 (18) 

e o e o 

For the example of Coulomb friction, the analysis can be simplified by expressing 
the, averaged properties in terms of several nondimensional parameters: (normalized 

quantities are denoted by a * superscript) 

k (B ) « k . T(B*) 09) 

co d e o 


c (B ) 

e o 




( 20 ) 


where 

c * 2m» * 2(m k ) 1/2 

o o d 


( 21 ) 



( 22 ) 


and 


B B k 

B * * -2- ■ -2-4 (23) 

ox.. aN 
crit r 


k tB *) and c *tB *) are the nondimensional averaged stiffness and damping of the 

CO CO 

friction element. These quantities are of particular interest because their values may 
be calculated in terms of the single nondimensional parameter, B q , the initial 
displacement divided by the displacement required for slip to occur. 
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For Coulomb friction from (4), (17), and (19), k * is given by 


k *(B *) 

e o 


1 

«r(B * - 1) 
o 



0(B*) - 0.5 sin[20(B*)3 



(24) 


where 


* 2 

cos0(B ) * 1 - -i , (B £ 1) 

B 

and need be calculated only once. The nondimensional frequency is then obtained 
from (12), (17), and (19) as 

=* ®.*> * ( r * V <0 )" 2 '«> 

d 

Lastly, from (5), (16), and (20), the average nondimensional damping is approximated 

as 


o 



which, upon integration, yields 


(26) 








♦ £n (B o # ) - 1 } 


(27) 


Note that the parameters are averaged only over the range in which they exhibit 
nonlinear behavior. The system is linear for B * less than one and the response can 

O 

be calculated using standard methods in that regime. 


The AA Method gives a linear estimate of the nonlinear behavior of the system. 
Results from the AA linearized systems are shown in Figures 2a-d as exponential 
decay envelopes, and are plotted with the Half Cycle Method and Runge Kutta results 
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for comparison. As with the Half-Cycle Method, we again see conservative 
estimates of system amplitudes and again the estimates are more precise for weakly 
nonlinear systems. 

An important result of the AA Method is shown in Figure 3. This figure shows the 

normalized quantities k *(B *), c *{B *), k *(B *), and c 1[B *) in terms of the 

eoeoeo eo 

nondimensional initial displacement B q *. It is observed that there are optimal system 

configurations which maximize either instantaneous or average damping in the system. 

These optimality conditions depend only on the nondimensional initial displacement 

B *. Recall from (23) that B * is dependent on the system parameters. Thus, by 
o o 

adjusting just one of these parameters, say, normal load for example, it is possible 
to optimally damp an existing system. For example, the points A, B, and C on the 
average damping curve can be thought of as three systems which are identical except 
for the tightness of the friction joint as indicated by the normal load (normal load « 
N, N A > N > N ). Comparing the linear average damping terms provided by the AA 
Method for the three systems, it is seen that system B, with normal load N g , has the 
highest value of average damping and is an optimally damped system for this set of 
basic parameters. The optimality of system B in an average sense was confirmed by 
fourth order Runge-Kutta simulations. The results for the three systems A, B, and C 
are shown in Figure 4, where successive extrema have been connected to form an 
envelope of decay. The normal load N_ is not simply the largest or smallest normal 

C 

load which could be applied, nor is it equal to the load which would be required to 
optimally damp only the first oscillatory swing In the transient motion (this normal 
load would correspond to the conditions of system D of Figure 3 which maximizes 
instantaneous damping). 

Runge-Kutta simulations confirm the optimality results which were readily provided 
by the AA Method. We observe, however, that the optimality result would not have 
been nearly so obvious if an exact time integration alone had been used to 
investigate this transient behavior. It is also noted that the AA Method provides 

general results in that it allows for parameter-based comparisons of different 
systems. Numerical time integration methods lack this generality, yielding instead 
results which are case-specific and thus more difficult to interpret when comparing 
systems. 

The AA Method leads to systems which are optimized, in an averaged sense, over 
the entire range of nonlinear behavior. The Half-Cycle Method may be used to 
generate systems which are optimal in a "first swing" sense. It is also possible to 
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optimize the system over other select ranges of nonlinear behavior. This is done by 
averaging the c ff and k equations over the particular range of concern, resulting in 
a new curve for c 1[B *). The optimal normal load (or other parameter) is the one 
which adjusts the nondimensional initial condition so as to maximize the value of 
this new c e ^B o *). The linearizations described so far make it possible to optimize 
the nonlinear system, in an average sense, over any range of nonlinear behavior. 

The Amplitude Averaging Method has yielded a general result which was not 
obvious from numerical time integrations. In the next section, the application of the 
AA Method to multiple mass systems is developed and it is shown that the 
principles of generality-of-results and optimal damping still apply, only that they now 
apply in a modal sense. 

Multiple Degree of Freedom Systems 

The transient analysis of multiple-degree-of-freedom (MDF) nonlinear systems using 
the AA Method can be accomplished by representing the system in modal form in 
terms of a sum of SDF nonlinear systems. The SDF components are linearized 
separately and are then combined to form a linear representation of the MDF system. 

In the linearization process, a nonlinear friction damper will be replaced by linear 
elements which approximate its behavior, as depicted in Figure 5. The original 2DF 
system incorporates a friction damper as a nonlinear element and will serve to 
illustrate the application of the Amplitude Averaging Method to MDF nonlinear 
systems. The equations of motion for the system may be written in matrix form as 

Mx + Cx ♦ Kx ■ - f 

n 

where 
x * 


M * 
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** 1 
k 2 + k 3 J 

MDF Linearization 

' i 

In the SDF analysis we were able to regard the initial amplitude across the friction 
damper as the maximum distance that the joint would be stretched during transient 
oscillation, in the MDF system, the initial amplitude across the joint is no longer 
guaranteed to be the maximum span which the joint experiences. Consequently, the 
nonlinear system cannot be linearized by simply considering the isolated friction 
damper and how much it is initially displaced. Instead, one must use a modal 
approach in the linearization process. 

It is necessary to decouple the nonlinear MDF system of Figure 5 into two SDF 
nonlinear systems so that the AA Method may be applied to each of the decoupled 
systems separately. Decoupling of linear systems is done routinely (see, for 
example, Thomson [7]) while nonlinear systems are not generally amenable to such 
analysis. The nonlinear decoupling and modal linearizations may be accomplished 
with the iterative approach summarized in Appendix A. In this approach, a set of 
converged eigenvalues and eigenvectors is obtained and used to form a modal 
(decoupled) representation of the nonlinear system. The decoupled systems are 
linearized using the AA Method and are transformed back to the original coordinates 
to yield a linear MDF system. 

Several approximations are made during the decoupling process. Viscous damping 
is generally small and is, therefore, neglected in order to avoid the inconvenience of 
dealing with complex eigenvectors. The complex formulation could be pursued in 
order to increase the overall accuracy of the linearization and estimation scheme. 
Another assumption made in the analysis is that each modal friction joint operates 
independently. This approximation significantly simplifies the calculations, but 
introduces an additional source of error in those cases where modes interact. This 
latter assumption is believed to be responsible for the fact that the modal AA 
response estimates are no longer conservative in some instances. 
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However, in general the results obtained from the AA Method compare reasonably 
well with results from direct time integration. The modal comparisons for two cases 
are shown in Figures 6 and 7. where the time integration solutions have been 
transformed from the original x coordinates into the converged modal coordinates, y 1 
and y 2 , by using the converged eigenvectors. Overall, the exact modal solutions 
{represented by the continuous curves) are seen to be similar in form to decaying 
sinusoids centered about a zero equilibrium state. The neatness (symmetry and 
sinusoidal appearance) of the exact results indicate that the converged linear 
decoupling does in fact represent the nonlinear behavior fairly well, as the response 
in x coordinates would transform poorly if the modal representation were not a 
reasonable estimate of system response. 

Several characteristics of the AA Method as applied to MDF systems are illustrated 
in Figures 6, 7, and 8. Firstly, the AA estimates are not conservative in the lowest 
mode when both modes are actively slipping. 2 Secondly, it was found that the larger 
the viscous damping and nonlinearities in a mode, the less "neat" the modal 
response looks (this is apparent in Figures 7b and 8b). Thirdly, a MDF system can be 
modally optimized in a fashion similar to that used to optimize the SDF system. As 
in the SDF case, the correct manipulation of the normal load shifts the 
nondimensionai modal initial condition to a condition which generates the optimal 
(modal) damping. It is noted that it is not generally possible to optimally damp both 
modes simultaneously, as the optimization of one mode results in a detrimental or 
non-optimizing shift of the other mode. Figure 6 shows a system which is not 
optimally damped in either mode. Figure 7 shows the same system after the normal 
load has been adjusted to optimally damp the higher frequency mode (mode 2). 
Notice that the results are plotted on different scales in Figures 6 and 7. In Figure 
8b, the extrema in the numerical time integration solutions for the non-optimized- 
mode 2 and optimized-mode 2 systems have been scaled for direct comparison 
(extrema plotted to form decay envelope). From Figure 8b it is obvious that the 
optimization has a significant effect on modal response. The same system was 
optimized in mode 1 and the results depicted in Figure 8a. Note that in this particular 
system the friction damping in the first mode was quite small compared to viscous 
damping and, consequently, optimizing joint damping had little impact on first mode 
response (an expanded scale was used to show the optimization result more clearly). 
The AA optimization results for both modes 1 and 2 were confirmed by Runge-Kutta 


2 

tt was found that if only one mode was actively slipping the accuracy of the AA method for MDF systems 
was comparable to that found for SDF systems. 
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simulations. Thus, although the modal AA extrema may be nonconservative in some 
cases, the AA indication of modal optimality remains accurate in the cases 

considered. 

The AA Method may thus be used to optimally damp select modes of particular 
concern. This is helpful in the design process in that it gives an easy indication of 
which systems are optimal in a given situation. In practice it may not be possible 

to obtain joints with characteristics that exactly correspond to this modeled 

optimality, in which case the issue may become one of selecting the most nearly 
optimal joint configuration from a variety of available designs. 

Conclusion 

This paper has discussed the specific application of the Amplitude Averaging 
Method to friction damped systems. However, the AA Method is a general 

linearization method applicable not only to friction damped systems but also to other 
systems which exhibit nonlinear hysteretic behavior. 

The AA Method is efficient and can be used to easily establish optimization 
conditions, subject to time integration verification. Familiar modal analysis may be 
applied to MDF nonlinear systems and systems may be optimized over specific 

ranges of nonlinear oscillation. Furthermore, the AA Method can be used to modally 
optimize MDF systems in order to suppress system response over specific frequency 
ranges. 

The AA Method may also be used as a comparison tool in the system design 
process. The physical parameters of the system may not be adjustable to the 
indicated optimal values. For example, in the case of jointed structures and friction 
damping, it is not generally possible to select physical joints with adjustable 
(optimizable) normal loads. In this situation the task may be one of selecting joints 
from a variety of designs. The AA Method yields the relative averaged damping in 
these designs for amplitudes of response which are representative of those 
encountered in practice and thus may be employed as a method of comparison in 
order to help choose the most nearly optimal design. 

The AA Method is an efficient design tool for two reasons. One is that the 
method is computationally efficient. In the cases considered in this paper, AA 
solutions could be calculated an order of magnitude more quickly than numerical time 
integration solutions. In addition, the method provides a system representation in 
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terms of linearized modes and, consequently, it becomes relatively easy to establish 
optimum system response. Again, it is our view that numerical time integration and 
the AA Method are complimentary approaches and that both have their places in the 
design of nonlinear systems. The AA Method is a computationally efficient approach 
which supplies analytical insight at the expense of accuracy, while the 
computationally intensive time integration approach provides verification and fine 
tuning of the results for select cases of interest. 
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A. ITERATIVE METHOD FOR DECOUPLING AND 
LINEARIZING A NONLINEAR MDF SYSTEM 


The goal here is to represent the friction joint by equivalent linear elements for a 
given set of initial displacements. To this end we must find a modal representation 
of the friction element which depends (as in the case of the SDF system) on the 
initial modal displacements. It is not possible, however, to calculate the mode 
shapes or eigenvectors of the system and find the initial modal displacements unless 
we know all the stiffnesses of the system (including those of the friction joint). 
This appendix summarizes an iterative procedure for simultaneously establishing the 
eigenvectors and initial modal displacements of the system and the equivalent joint 
properties of the friction element. 


This method neglects damping in determining the eigenvectors of the system. This 
is a reasonable approximation for lightly damped systems. There are seven steps in 
the process, which are illustrated by applying them to the two-degrees-of-freedom 
(2DF) system of Figure 5. 

1. Write the governing equation of the system in matrix form, neglecting 
viscous damping, e.g., for the 2DF system of Figure 5 


Mx + K x = f 
o n 



(A.1) 


where K * K. 
o 

2. Find the eigenvalues and eigenvectors of the linear part of the system, 
momentarily disregarding the nonlinear term f . Form the matrix of 

eigenvectors, U, e.g., n 

K u . * X.Mu. eigenvalues X.,X„ 

o i i i 3 12 

eigenvectors u n * 

U * [u ,u 2 ] ... matrix of eigenvectors 


-t] 


(A.2) 


3. Assume x can be represented in terms of modal coordinates y * (y^y^, 
i.e.. 


x > Uy 


(A.3) 




then 
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y = U \ 

The modal initial conditions can be found using (A.3) and the orthogonality 
of the eigenvectors (with respect to M) 

y = U _1 x 
o o 

4. Using x = Uy in (A.1) yields 

MUy + K Uy = f 
o n 

U T MUy + U T K o Uy = U T f 


-1 

1 


5. Find linearized stiffnesses for each independent mode using the AA 
method. This step yields 

U T MUy ♦ U t K Uy * - K y - C y 
o' q' q 

(Analogous to: 


r r . 

-f « — -x + — x^>-kx-cx 
n B wB e e 

as done in the SDF case in equations (7,8,16,17,18) in the main body of 
the paper.) K is a diagonal matrix, each diagonal element is the 
linearized stiffness of the friction element in its respective mode, e.g., in 
the example 2DF system, slip occurs when 

,N . k d (x, - X 2 > ■ k a y, * k dj y 2 

where 

k d - k d u. [1,-13 

i 

Assume that slip occurs independently in each mode. Then the modal 
displacement required for slip is 


//N /<N 

Yj * — (analogous to x crjt « — ) 
crit 
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The modal initial conditions were found in Step 3. Given the initial 
conditions and the slip conditions, the AA method may be applied to each 
mode separately in order to linearize the decoupled systems. From the 
first mode, the linearized stiffness term is -k . Mode 2 yields -k . 

e i ®2 

fk 0 1 

Form K : K * e l 

Q Q L 0 k J 

e 2 J 

The linearized modal stiffnesses have been found (C is found 
simultaneously but not used at this point in the analysis). 

6. Convert back to original coordinates: 

ii _1 

x * U y 


Substituting into the linearized governing equation from Step 5: 


-1 


Mx 


K x 

o 


- U' 


K U -1 x 


(neglect C q ) 


define 


K 


mod 




then 

* * + “mod* ■ 0 


7. 


K ♦ K 


mod' 


then 


Mx ♦ KjX ■ 0 


But this system has different eigenvalues and eigenvectors than the 
original system considered in Step 2. The linearization of Step 5 was 
dependent on the eigenvectors of the system, which have changed after 
linearization. To obtain a more accurate linearization, repeat Steps 1-7 
using the eigenvectors from the M 1 ,K 1 system. 


The process (Steps 1-7) may be repeated as many times as needed to obtain a 
converged set of eigenvectors using M,K (j) . In practice, this method usually 
converged in four iterations and never required more than six iterations. 


While the linearized damping can be readily calculated at any iteration, it is 
neglected along with any original viscous damping in estimating the eigenvectors of 
the system. For this reason, the method is less accurate when high viscous damping 
or large nonlinearities are present in the system. A complex eigenvalue and 
eigenvector approach could be used to improve the accuracy when the damping terms 
are large. 
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In the current approach the equivalent linear damping from the friction joint is 
calculated for each mode only after the eigenvectors have fully converged. This is 
added to the modal viscous damping to establish the total effective damping 
(analogous to (16)) for each mode. The response of each mode is then calculated in 
the same manner as for a linear single degree of freedom system to obtain the 
comparisons indicated in Figures 6 and 7. The physical displacements may be readily 
calculated from (A.3). 



20 


Figure 1 
Figure 2 
Figure 3 
Figure 4 
Figure 5 
Figure 6 

Figure 7 

Figure 8 


Figure Captions 

Single Degree of Freedom System. 

Comparison of Results for m * k + k = uN/k . * 1 and c * 0.02. 
Instantaneous and Average Values of Stiffness and Damping. 
Response Envelopes Indicate Optimum Damping. 

Schematic of Linearization Process. 

Modal Response of MDF System: k * k * m * 1, k * 2, c * c 
* 0.01, k d * 0.25, and /iN * 0.1. 

Modal Response When Second Mode Optimally Damped: fiN 

Changed to 0.35. 

Comparison of Modal Response Envelopes. 
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1. INTRODUCTION 

Many engineering structures are subjected to motions in which limited slip can occur 
between members. This is frequently the case for bolted and riveted structures and can also 
occur in joints found in space applications. Truss structures designed for space often are 
constructed with joints that can rotate and lock into place. This allows the structure to collapse 
into a compact package that is easily transported to space and then expanded. Typically, this 
type of joint mechanism is subject to limited slip behavior because of "slop” due to machining 
tolerances. In addition, joint slip is usually constrained by friction forces within the mechanism. 
This means that the joint’s motion is controlled by friction until it reaches the end of the free 
play and then it is constrained elastically. This type of nonlinear behavior can have a strong 
influence on system response. In particular, because of the limited slip or the ’gap’ aspect of the 
joint behavior scaled models of the jointed structure may not respond dynamically in the same 
manner as the full scale structure in space. This may be the case for two reasons. 

The first reason that a scaled model of a structure may not have the same dynamic 
response is that tolerances in the joints are not scaled. Typically, engineers wish to have joints 
that behave as linearly as possible. Consequently, when they are designing the full-scale joint 
they stipulate as small tolerances as possible to reduce slop and the associated nonlinear behavior. 
When sub-scale models of the joint are developed for testing, the same absolute tolerances are 
stipulated since they cannot be further reduced. As a result, the free-play in the sub-scale joints 
is, on average, significantly greater than in full-scale joints and the model structure will exhibit 
stronger nonlinear behavior. 

A second reason that a scaled model of a structure may behave differently is that it may 
be tested under gravitational loading that does not occur in space. While there are laboratory 
procedures for supporting two dimensional structures ( i.e. structures that have modes in which 
the motions lie in a plane) in such a way as to negate gravitational effects, they are not readily 
applicable to more complex, three dimensional structures. As a result, it is in general not 
possible to eliminate gravitationally induced preloads in the joint. The preload introduces a static 
displacement in the joint that tends to eliminate the free-play in one direction. Thus, for 
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example, the joint may behave elastically in the direction of the applied preload, and slip only 
when motion occurs in the opposing direction (in contrast to a joint without any static loads that 
slips in both directions). As a result, from a dynamics point of view, a joint with preloads 
tends to be effectively stiffer and contributes less damping to the system. A goal of the current 

research is to gain a better understanding of this effect. 

This research considers the dynamic response of the lumped parameter system depicted in 
Figure 1. The system may be viewed as either a model of a single degree of freedom oscillator 
or a single mode representation of a general structural system as characterized by a given modal 
mass, stiffness and damping. The nonlinear element that represents the joint is indicated by the 
linkage having the spring stiffness k . In this linkage the friction element can slide once the 
magnitude of the force builds up to fiN. Slip at the friction contact point is limited by the stops 

located at distances b ] and and the total gap is given by <5. The system is subjected to a 

sinusoidal excitation and a static preload f . This system was chosen for examination because it 
is the simplest system which exhibits the characteristics of interest and yet, in a modal sense, has 
implications to a broad class of systems. To date we have concentrated on understanding the 
steady-state response of this system. 

The steady-state response of the system is important for two reasons. The first is that it is 
easier to analyze than the transient response and, consequently, it is easier to establish important 
nondimensionalized system parameters. Since these parameters control transient response as well 
as steady-state response the final issue is how much they change in going from a full-scale to a 

sub-scale system. For example, one suspects that static preloads have little effect when the 

system’s response is large since the amount of energy dissipated by friction is then independent of 
the static loads. This can be readily quantified by analyzing the system’s steady-state response to 
establish for what range of dimensionless parameters this is true. This analysis is easier to 
perform for steady-state response because well established techniques (the describing function 
approach) exist for determining approximate solutions. The resulting analytical solutions are 
superior to numerical simulations because they are easier to compute and more general. Secondly, 
it has been shown in other work recently performed at Carnegie Mellon that the steady-state 
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behavior of a system can be used to estimate its transient response. Thus, it is clear that gaining 
an understanding of the factors that control the steady state response also has a direct bearing on 
how it will behave in the transient regime. 

I wan has analyzed a similar limited slip system in [1]. It differs from that considered 
here in that it did not contain the spring with stiffness k (refer to Figure 1) nor did he consider 
static preloads. Both of these factors are important in space applications. 

2. ANALYSIS 

The solution approach is explained in this section. In Part 1, the governing equation for the 
system is described. Its normalized form is given to provide a basis for parametric studies. The 
approximate method is introduced in Part 2, and the equations for calculating the frequency 
response are derived. Part 3 presents the procedures for obtaining the resonant responses. Part 4 
shows how the frequencies and amplitudes where the vertical jumps occur are obtained. The 
procedures described in these sections provide an approach for understanding the unstable response 
phenomenon discussed in the latter sections. 

2.1. Formulation of Governing Equations 

The equation of motion for the system shown in Figure 1 is, 

mx(t)+cx(t) + k x(t) = 1 cos wt + f - f (x,5,,5J (1) 

o I n 1 2 

where /, is the static load, and f is the nonlinear friction force at the joint. 

1 n 

The contact pressure on the joint, N , is assumed constant, and the inertia of the friction 
element is neglected. Therefore, f n is constant throughout the length of the friction element. 

The friction element will slide once the magnitude of the nonlinear force f equals N 
multiplied by /i, the friction coefficient of the materials in contact. 

/iN 

The nonlinear force f is a function of the mass displacement X. When | x | £ — , where 

k d 

k is the stiffness of the friction element, the system is linear. When | x | £ — , and the stop 

k d 

distances 3 , c> 2 are sufficiently large, the system has elasto-plastic behavior. When the stop 
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distances, <5 ] , <5 , are not large enough, the system experiences a sudden hardening phenomenon. 
If the mass displacement x is periodic, several different hysteresis cycles of 1 ^ vs. x can be 
drawn, as shown in Figures 2-6 ( In the figures, A is the amplitude of the steady-state vibration 
and B is the static offset). 

To render a foundation for the parametric study, equation 1 can be converted into a 
nondimensional form as shown below 

d 2 x/ dt 2 + C dx/ dt + x = cos wt + \ ] - X 2 f( x, S 2 /o } ) (2) 

Where, 

x = x/x and x = f /k 

t = t/t o and t o = x/m/k 

X, = ljf n , t = 

X, = pN/ f o , w-w/t o 

As can be seen, the nondimensional scales are unities if the quantities m , f were set to 

o 

unity. This fact is adopted in the sections below as a way to present general results. 

One purpose of this study is to understand the steady-state response of the system shown in 
Figure 1. Although this task can be accomplished by solving equation 1 numerically for a 
sufficiently long period of time, such an approach is inefficient and computationally expensive. 
Therefore, an alternate, more efficient method is developed. The direct, long time solution is 
then used only to verify the accuracy of the analytical results in certain representative cases. 

2.2. Approximate Method 

Based on the fact that the excitation is periodic; we assume that the response is also 
periodic. Therefore, the displacement x can be expanded by Fourier series. In the steady-state, the 
response is approximately harmonic, i.e. 

x - Ac os d + B , where d = wt - f (3) 

In the above, the offset term B is included to compensate for the shift due to the static load / . 
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The method used here was adopted in several earlier papers [1] [2] [3] that studied the 
friction damping of structures. Since the nonlinear force 1 depends on the oscillation of 
displacement x, it is also harmonic and can be linearized by representing it in terms of a 
truncated Fourier series as shown below. 

/=/.+/ cos 6 + f sin 6 (4) 

n b C S 

Where 

?NFJa,p,y.q) 

f c - pN F{ a , /S , y , q ) (5) 

~ F ( a , P , y , q > 

In the above, F , F , F^ are nondimensional Fourier coefficients established in the Appendix and 
a, y and q are dimensionless parameters also established in the Appendix (see equation 43). 

It is noted that since the dynamic displacement is proportional to cos 6, it is apparent from 
equation 4 that contributes to the dynamic stiffness of the system while provides damping. 
The stops shown in Figure 1 suddenly contribute to the elements stiffness and significantly 

inhibits joint damping. As can be seen in Figure 7, the magnitude of f decreases with the 

amplitude right after the limit stops are hit. 

To obtain the approximate solution, equations 3 and 5 are substituted into equation 1. 

After equating the linearly independent terms, the phase shift f is eliminated by using the relation 
sin 2 f + cos 2 f * 1. Then the following two nonlinear algebraic equations are obtained. 

kB+f h -f=0 = g,(A,B,w) (6) 

Dl 1 

Uk-mu> 2 )A + f c ] 2 + [wcA - f s ] 2 - Pm 0 =g 2 ( A , B , w) (7) 

Equations 6 and 7 are solved iteratively to get the frequency response of the vibratory 

amplitude A and the permanent offset B. 
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2.3. Calculation of Resonant Response 

By holding the system parameters constant, the resonant vibratory amplitude A is obtained 
as 


SA 

— =0 at w — w 

Sw 


( 8 ) 


To relate equation 8 with equations 6 and 7, it is observed that A and B are both 
functions of w. Consequently, g , g_> are only functions of w, and 


d l±- + o 

Ba Bw Bb Bw Bw 

d 9 2 3 9 2 BA 9 2 dB $9 2 

d w BA Bw BB Bw Bw 

From equation 6, we know Bg } /Bw = 0 , thus equations 8 and 9 imply 
SB 

= o and since from equation 6 * k # 0 

SB Sw SB 

We have 

Bb 

— = 0 at w » w 

Sw 

Consequently, from equation 10 


(9) 


( 10 ) 


( 11 ) 


( 12 ) 


Sw 


- {clwcA- f s ]-2mw[(k- w 2 m)A + f c ’]}(2A)-0 at w - w n 


(13) 


Thus, the three equations that determine the resonant values A and B at the corresponding 
resonant frequency of excitation are 
9.( A. B. w) = 0 

l m m m 

gj A.B.w m )*0 (14) 

Bg^ 

— ( A m , B m , w m ) = 0 

* m m m 

Ow 

To solve equation 14, it is noted that equation 6 is decoupled from equations 7 and 13. 
And from the relation (7) + (13)«; = 0, we have 
Fw 4 + Gw" + H — 0 

m m 


( 15 ) 
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where 


F = 5( mA )' 

m 

6=3 {( C- ~ 2mk )A 2 -2mA f } 

m me 

H = (A A ) 2 + f~ + 2k A f +f 2 -r 


(16) 


me s o 


By substituting equation 13 into equation 7, we have 


{w cA -f } 2 { 1 + ( > 2 } = f 2 

2 mw 


mm s 


(17) 


Therefore* from equations 14, 15, 16 and 17, the following iterative equations for 

calculating the resonant response can be obtained 

B 

m k 


A m s{f S\ f o\ / ) / cw a 


(18) 


1 


<* = {(-)(- 6+ /6"-4/7/ )} m 

m 2p 

To solve equation 18, group values of ( A , B , w ) are first initialized to calculate the 

^ m m m 

Fourier coefficients /, 1 , i . Iterations are then made utilizing equation 18 until A , B , and 

b s c mm 

converge. 


2.4. Calculation of Vertical Jumps 

When the system parameters are constant, the vertical jumps occur when 

BA Bw 

— -400 or — = 0 at w = w 


Bw 


Ba 


(19) 


The above constraint is imposed on equation 7 to get the following relation as 

dg 2 ^9 2 Bw d9 2 

= + — s — = 0 at w-w 

dA $a Ba Ba 1 


where 


« *1 h 


Bf_ 


= [ (k-mw 2 ) A + f ]( k-mw 2 + — ) 

Ba c Ba 

+ ( - f )( wc ) = 0 at w-w 

5 Ba j 


( 20 ) 
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From the relation ( ABf /BA + >< 20) - (2Af)(l) - 0 , we obtain the following equation 

Rw A + Sw~ + T = 0 (21) 


where 

3f. 

= ( A — 1 - f )(mA) 2 
BA s 

Bf f 3/^ 3f s 

S = { (c--2mk)( A— - f ) + (2m)(f--f — )}A 2 (22) 

s s 3>* c a<4 

3f 3/ c a/ a/, 

T= {(Ak -If )—-<* + 2 — ) / }(Ar <4 2 ) - (2/tf ></ — + / — ) 

3/4 3/4 '3/4 3/4 

3/ 

+ (f 2 + f 2 -f 2 )(f + A — ) 

S C O S - 


and the frequencies at which jumps occur are obtained by 

{( -S±v/S"- 4/?r )/ 2 R)' n (23) 


Analogous to the approach used in the last section, equation 23 together with equations 6 
and 7 provide the basis for calculation of the vertical jumps. 


3. NUMERICAL RESULTS 

The numerical results are presented in this section. Part 1 shows representative plots of the 
amplitude as a function of frequency of excitation and discusses how various parameters affect 
the response. Part 2 indicates how changes in the system parameters affect the peak response 
amplitude. The effects of gap length on dynamic response are shown in Part 3 and Part 4. Part 
3 describes the existance of disconnected amplitude solutions and how .the gap and static preload 
affect them. Part 4 introduces the minimum gap length for unstable response to appear and how 
instability is affected by preload. Finally, two cases of transient response are given in Part 4. 
They show that the system parameters that dictate the steady-state behavior of the system also 
control transient response. 
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3.1. Frequency Response 

To examine the validity of the approximate method; equation (1) is solved directly using 
the Runge-Kutta method. The steady-state response is then obtained from the "long time" transient 
solution. The results are compared with the ones calculated by the approximate method. These 
comparisons are plotted in Figures 8-1 2 in which the discrete symbols denote the transient 
solutions. These figures show that a good agreement exists between the two approaches. 

In Figure 8, the influence of gap length <3 on the vibratory amplitude frequency response is 
shown. For <5=0, the system is linear. For 5 oo, the friction joint is dominated by friction 
damping, and the system behaves plastically. In between, a transition region exists for a certain 
range of <5 values where the system shifts from elastic to plastic behavior. In this region, the 
vibratory amplitude A may become multi-valued at some fixed excitation frequencies. This 
phenomenon is typically referred to as "unstable response". A typical case for <5 = 2 0 is shown in 
the figure. 

Unstable response can also be seen in Figure 9, in which the amplitude is calculated for 
several different values of the normal load ///V. For pN = 2.0, the system response is stable 
(vibratory amplitude is single valued). Unstable behavior is clearly seen in cases of pN = 1.0 and 
1.5. For pN = 1.5, the vibratory amplitude has two disconnected response curves. For pN = 1.0, 
the response curve becomes connected, even though point jumps remain at frequencies between 
1.035 and 1.04. This implies that for even lower pN values, the system response will become 
stable again. The physical explanation of this behavior is that at low pN values, the friction 
damping has little effect on the system response which is primarily dominated by the constant 
stiffness k and viscous damping c, and consequently, the system behaves linearly. For high pN 
values, the friction element is partially stuck due to the stronger friction resistance. Under this 
condition, neither of the limit stops is hit during the motion cycle. As a result the energy 
dissipation per cycle monotonously increases with amplitude, and the system response is stable. For 
normal loads in the median range, as the cases of pN = 1.5 and 1.0 shown in the figure, the 
limit stops are hit by the element during the motion cycle. The system response becomes meta- 
stable or unstable. This is due to the fact that the energy dissipated by friction per cycle 
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remains constant while, because of the limit stops, the element stiffness significantly increases 
with amplitude. This results in the multi-valued response indicated in Figure 9. 

Multi-valued behavior is also seen in Figure 10; where the amplitude is calculated for 
various stiffness ratios k & / k. A region of unstable behavior is seen for k / k = 0.25, where a 
disconnected closed curve occurs. The appearance of disconnected regions in the more general 
system model analyzed here agrees well with those observed in Iwan’s [1] work. 

The effect of static load f on the vibratory amplitude frequency response is shown in 
Figures 11 and 12. The static load causes a permanent offset represented by the symbol B, as 
indicated in equation 3. When the system is excited by a large excitation force f at the 
resonant frequency, both limit stops ( shown in Figure 1 ) are hit during the motion cycle, and 
the offset has little effect on the vibratory amplitude. However, at some off-resonant frequencies, 
the response is sufficiently low that there is a difference in the response since without the static 
preload the limit stops would not be encountered during the motion cycle. However, once the 
static load is applied, the permanent offset thus generated adds to the dynamic displacement and 
causes the friction element to hit one of the stops and reduces damping in the system. 
Consequently, the vibratory amplitude increases at off resonant frequencies as illustrated in the 
figure. 


Based on the above argument, it is then understandable that the static load may even raise 
the maximum amplitude of the system under low excitation force f . This fact is depicted in 
Figure 12 where the maximum amplitude is increased when the static load f = 10. 

Thus, it is seen that under certain conditions static preload can have a significant effect on 
vibratory response. For some excitation, the presence of a large static load may even eliminate 
the multi-valued response and stabilize the system. This means that gravitational effects may 
strongly affect the dynamic response of structures that are ground tested and it might be quite 
different from the actual response that would occur in space. Therefore, this effect needs to be 
carefully analyzed in order to properly interpret laboratory test data so as to predict the response 
of the full-scale structure in space. 
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3.2. Peak Amplitude Response 1 

First consider the system’s response when there are no stops to limit slip in the nonlinear 
element. The amplitude of peak response will then depend on how stiff the nonlinear spring k d 
(relative to the system’s stiffness, k) and the value of the slip load, pN. Representative results 
are presented in Figure 13. For sufficiently large /*/V, the friction element is stuck and the 
system is linear. The system also behaves linearly when pN is zero since no energy is dissipated 
by friction. It can be seen from the plot that the friction damping effect is only significant for 
///V between 0 and some maximum value, Furthermore, a minimum peak vibratory 

amplitude exists for a certain pN value, as shown in the figure. This fact has been discovered 
and experimentally verified by Griffin [3]. 

Figure 14 shows the plot of resonant frequency vs. normal load ///V at gap length <5 = oo. 
It is observed that resonant frequency increases with the normal load, and hits a maximum when 
the normal load is sufficiently large that the friction element is always stuck. Again this 
indicates that the effect of normal load is only significant over a finite range. 

The effect of gap length on the peak vibratory amplitude response is indicated by the 
results depicted in Figure 15. A plot of the resonant frequency vs. gap length is given in Figure 
16. In both figures, the y in the abscissa denotes the nondimensional gap length introduced in 

equation 43. In both plots it can be seen that for constant normal load // N , a transition region 

exists in which the gap length strongly affects the dynamic response. For gaps larger than a 
certain value, the limit stops shown in Figure 1 will no longer limit slip in the friction element 

and, consequently, further increases in gap length will not affect the peak response. It is clear 

that the gap length effect is only significant in the transition region. This is depicted in Figures 

17 and 18, which are the scaled versions of the Figures 15 and 16. The critical gap length c> cr 

is defined as the value beyond which the peak response is unchanged. This value which will be 

investigated in the next section is used as the scale length in the abscissa for both plots. 


Since in an unstable response condition the system may have several resonant amplitudes, here to avoid ambiguity, "peak amplitude** 
refers to the maximum resonant amplitude. 
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At S-0 the system peak response is dominated by stiffness and viscous damping. At <$=oo, 
the system peak response depends on the stiffness, and both viscous and friction damping. The 
difference between the two amplitudes and frequencies are used to scale the ordinates in Figures 
1 7 and 1 8 respectively. 

Since the peak response is only affected by the normal load pN = 0 to ///V , the results 

for different normal loads all lie fairly closely together in a band as indicated in Figure 17 for 
amplitude and Figure 18 for frequency. It was shown in Figure 13, that the lowest peak 
response comes at a specific pN value (0.1 for the system considered here) and this case 
provides a lower bound on the response when the data is depicted as in Figure 17. The normal 
loads which are greater or less than this value then approach this common lower bound of the 
band. 


As the normal load goes to zero or when it approaches /'A/ max the system becomes linear 
and the stiffness has no effect on the lower bound curve of the band. However, an increase in 
the nonlinear element’s stiffness kjk can affect the upper bound of the curves. This is shown in 
Figures 19 and 20 for peak amplitude and frequency respectively. 

3.3. The Effect of Critical Gap Length on Unstable Response 

The phenomenon of disconnected amplitude solutions is examined in this section by 
considering the equations for maximum response. First, equation 17 is rearranged into the 
following form. 

-/ < A,B)* \fj /Vl- + (c/lmw'j* - w m cA * g{ A,B ,w m ) (24) 

The offset B can be calculated from the first equation in 18 as a function of A and then 
w m from the last equation in 18. Therefore, conceptually, functions g s and f are dependent on 
the amplitude A only. Thus, these functions can be presented in a plot with amplitude A as the 
common abscissa. These are shown in Figures 21 and 22 for a range of values of f and f } 
respectively. In the plots, the intersections indicate the resonant amplitudes. 


As can be seen, the function g ! is almost a linear function of A due to the fact that the 
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resonant frequency is nearly constant ( approximately 1.0 ). In Figure 21, it is clear that for 

different excitation levels, an unstable region of two disconnected amplitude response curves exists 

between points cl and c2. This is indicated by the three resonant amplitudes shown in the 

figure: the upper two values are the maximum and minimum amplitudes of the upper closed 

response curve, while the lowest one is the maximum amplitude of the lower response curve. 
Physically, point cl indicates the maximum excitation for which the stable response is obtained. 

In this situation, both of the limit stops shown in Figure 1 are not encountered by the friction 

element during the motion cycle. This indicates that for constant excitation a critical gap 
length <5 cr exists beyond which the response remains unchanged. This critical value can be 

obtained by the geometric requirement at point cl that the slopes of the curves as well as their 

values are equal, i.e. 

— = wc and ~f = g , when B is constant (25) 

BA BA 


Point c2 indicates the maximum excitation for disconnected solutions to appear. 

Alternatively, for constant excitation f ^ this point indicates the maximum gap length 5^ for 
which the system has disconnected amplitude solutions. It is clear from the figure that in the 

disconnected amplitude response region, an upper response curve exists ( indicated by the upper 
two resonant amplitudes ) when the friction element hits the stops during the motion cycle, while 

a lower response curve exists ( indicated by the lower resonant amplitude ) when it does not 

reach at least one of the stops. Since the lower response curve is not affected by the gap ( 
stops are not hit ), the lower resonant amplitude is the same for gap lengths larger than or equal 
to 5 at a constant excitation f Q . As can be seen in the figure, point c2 is where the friction 
element touches both of the stops at the lower resonant amplitude. Therefore, the gap length <5 m 


can be determined from the equation below 


5 

m 


= 2 (A 


oo 



where A^ - the resonant amplitude at 5 =00 


(26) 


In the above, the A is used since the lower resonant amplitude is the same for gap length 
larger than or equal to i5 fi) . 
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The effect of static preload on the region where the disconnected solutions occur is 
depicted in Figure 22. It is apparent that the static load tends to play an on/off role on this 
region. When this load is below a certain value, only the lowest resonant amplitude is affected 
and the response remains disconnected. Once this load is raised above a critical value ( 20 in 
the case showm ), the amplitude solutions become single valued. 

3.4. Minimum Gap Length for Vertical Jumps 

As was discussed in the last section, disconnected solutions exist between limit gaps t> cr and 
<5 m for constant excitation / . For gap length less than <3 m , the system response becomes 
connected; however, discrete point jumps remain. Thus, unstable response can occur until the gap 
length is further reduced below a limit <5.. 

The reponse depicted in Figure 23 may be used to illustrate this argument. The locus of 
vertical jumps is presented in broken lines which define a closed region. In this plot, the 
excitation f is varied while the gap is kept constant. The intersections of the response curves 
and the jump locus indicate where individual point jumps occur. In the figure, it is clear that 
vertical jumps only occur when the limit stops shown in Figure 1 are hit during the motion 
cycle. The jumps disappear when the excitation f is large enough that the response curve passes 
beyond point a ] or a 2 ( whichever is passed by secondly ). At constant excitation, the limit gap 
2. for vertical jump can be obtained when the response curve intersects the jump locus at either 
point a } or a 2 ( whichever leads to the smaller 5. ). The choice between points a } and a 2 as 
the passing point is affected by the stiffness * d / k. When the stiffness k d / k-> oo, the response 
curve is more slanted to the left and point a 2 will be passed by first in decreasing gap length at 
constant excitation. The point a } is then the intersection used to determine the limit gap and 
vice versa. 

In Figure 24, the effect of static preload on the jump locus is depicted. It is seen that 
the closed region of the jump locus is shrinked as the load is increased. At the same time, the 
jump limit gap S. is raised in value. 

The limit gaps 2 cr , <5 m , 2. as functions of the normal load pN are shown in Figure 25. 
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As shown in the plot, three regions exist. When the normal load pN approachs zero, regions II 
and III shrink, and all the limit gaps 3 cr , 3 m , 3. tend to infinity. The system response then 
becomes linear. 

At a large value of /*/V, limit gap curves 3 cr and <$ m superpose on each other, and the 
disconnected response ( region II) disapears. However, the discrete jumps (region III) remain. 
This remains the case until uN reaches uN , after which the friction element becomes stuck 

r r max 

and the gap no longer affects response. As can be implied from the figure, the limit gap curves 
3 cr , 5 m then descend to merge with d. at one point. 

3.5. The Effects on the Transient Solutions 

One reason for studying the system’s steady-state is to gain a better understanding of its 
transient response. We have been examining conditions under which the steady-state response is 
unstable. These conditions also affect transient response. For example, when the system 
parameters are such that its response tends to be unstable then slightly altering one of the 
parameters can lead to quite different transient response. 

Figures 26 and 27 indicate how the gap affects the transient response. Figure 26 shows 
the initial response of the system while Figure 27 shows how the vibratory amplitude changes 
over a large period of time. Two cases are plotted in each figure, one corresponding to 3-23 
and the other to 3 = 30. For short times ( Figure 26 ), responses are identical. However, for 
longer times they differ radically. It is known from the previous sections that for conditions 
corresponding to unstable response, different initial conditions may lead to different vibratory 
amplitude at the steady-state. In the previous section, we showed that the gap length also affects 
the unstable response. Consequently, even for the same initial conditions, changes in the gap 
length can lead to quite different transient responses. 

The effect of static preload is illustrated in Figures 28 and 29. Again, Figure 28 shows 
the plots of the displacement for short times, and Figure 29 shows the long ..time response. In 
this particular case, the load accentuates the response for short times, and attentuates it later. 
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In regard to the ground testing on a sub-scaled model of a large space structure, the above 
plots indicate an important result. The examples shown in Figures 27 and 29 imply that a larger 
gap may reduce the amplitude of response while the applied static load may contribute an 
opposite effect. Compared to the full-scale structure, the sub-scaled model tends to have larger 
joint tolerances relative to its smaller overall size. In addition, a large static preload ( from 
gravity) may be present during ground testing. Since these two factors sometimes contribute 

opposite effects, they have an unknown qualitative and quantitative effect on the test results. 
This indicates the need for a careful analysis of the experiment data in order to predict properly 
the response of the full-scale structure in space. 

4. CONCLUSION 

The approximate solution method based on a describing function approach has been verified 
to be an effective tool in studying the response of systems containing a friction joint with limited 
slip. It is believed that this approach could be extended to investigate the steady-state response 
of a large truss-like jointed structure and could provide a basis to calculate transient response. 

It has been found that under certain conditions the vibratory amplitude of the system may 

become multi-valued at certain excitation frequencies. A large jointed structure may have 

significant variations in the joint properties from joint to joint due to machining differences. As 

a result some of the joints may experience the conditions that lead to multi-value response. 
Consequently, it may be difficult to get repeatable experimental results since the structure may 
settle into different patterns of response depending on fairly subtle aspects of how the loads are 
applied. The difficulty is intensified for laboratory tests of sub-scaled models of space structures 
since the scatter in the joint’s free-play is relatively larger (with respect to the smaller 
dimensions of the model joint). 

This work also shows that in some circumstances, the presence of static preloads may affect 
significantly the system’s dynamic response. In this case if a large structure is ground tested it 
may be necessary to first correlate the data with an analytical model of the structure that 
includes preload effects. Then eliminate the preloads in the model and predict its behavior in 


space. 
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Lastly, a way of looking at system response in terms of non-dimensional lzed joint properties 
has been developed. These curves are useful because they indicate how sensitive a system is to 
changes in joint tolerances. When one makes a scale-model of a joint the relative tolerances and 
free-play in the joint increase. The computed results indicate an approach for assessing system 
sensitivity to such changes. 


5. FIGURES 
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Figure 1: Single Degree of Freedom System with ■ Limit Slip Joint 
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Figure 2: Non Slip 
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Figure 3: Pure Slip 




Figure S: Slip with Lower Limit 


Figure 6: Slip with both Limits 
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Figure 8: Vibratory Amplitude vs. Frequency with Varied Gap Length b, 

f <T 1 ' f x a °> ™ = hc =0.01, k d / k * 0.1,/iN * 1.0, q * 0.5 
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Figure 9: Vibratory Amplitude vs. Frequency with Varied Normal Load ft 

1.0, /,* 0, m */, c = 0.0 f, kjk “0.1, J«20, <7«0.5 
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Figure 10: Vibratory Amplitude vs. Frequency with Varied Stiffness kjk , 

^,*0. m = 7, c = 0.07, ii /V* 1.0, 3*20, q = 0.5 
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Figure 12: Static Load on Amplitude Frequency Response Under Low Excitation. 

f a - 0 . 3 , m= 1 , c*0.01,kjfk*Q.Ql, 1.0, a * 20 . q ~ 0.5 
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NORMAL LOAD. pN 

FlflUr * ,3: Max - Vibratory Amplitude vs. Normal Load pN when 

Gap Length i * oo, q • 0.5 
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Figure 14 : 


ResoMnt Frequency vs. Normal Load pN when 
Gap Length i » oo, g * 0.5 
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Figure IS: Max. Vibratory Amplitude vs. 
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Figure 16: Resonant Frequency vs. y, q - 0.5 






Scaled Max. Amplitude vs. <5/<5 cr with Varied Normal Load /iN, q =0.5 
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NONDIMENSION GAP LENGTH, l / b v 

Figure 18: Scaled Resonant Frequency vs. l/l u vith Varied Normal Load f N, q m 0.S 
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NONDIMENSION GAP LENGTH, b / b 

Cf 

Figure 19: Scaled Max Amplitude vt. b/S with Varied Stiffness Q 9 

*' k cr 
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VIBRATORY AMPLITUDE A 


Figure 22: The Effect of Static Preload on the Resonant Amplitude Response, 

f o *1.0, m *1.0, k =1.0, k* 0.1, c = 0.01, /i/V=1.2, 5*20, <7*0.5 
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The Limit Gap Curves as affected by the Normal Load pN, 

1 *1.0. /,= 0, ‘A:* 1.0, m=1.0, k *0.1, c*0.01, £«20, <7«0.5 
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Figure 25: 
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A 2*23.0 time 

— 2 * 30.0 

Figure 26: Effect of the Cap Length on the Transient Response [I] : Displacement vs. Time, 

initial pos. 1 0 , initial vel. * 30, excited at u m 1.04 

f*1.0, f,*0, Xr— 1.0, /n* 1.0, *j®0.1, c*0.01, /i/V* 1.2, <7*0.5 
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Figure 27: Effect of the Gap Length on the Transient Response [II] : Amplitude decay 

initial pos. =0, initial vel. * 50, excited at «=1.04 
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Figure 28: Effect of the Static Load f ^ on the transient response [I] : Displacement vs. Time, 

initial pos. = 0, initial vel. * 0, excited at w=1.06, 

/ o *=1.0, m= 1.0, *=1.0, *=0.1, c=0.01, ///V=2.5, 5=6.0, <7=0.5 
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Figure 29: Effect of the Static Load f x on the trasient response [II] : Amplitue Decay, 

initial pos. e 0, initial vel. s 0, excited at a?* 1.06 
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I. APPENDIX : Calculation of Fourier Coefficient of Nonlinear Force f 

II 

The variations of the nonlinear friction force 1 is considered as follows: 

n 

The friction contact sticks initially, and / varies linearly with slope k . When f = ± uN , 
it starts to slide , and f remains constant until the displacement x achieves an extremum or the 
friction contact hits either one of the stops. In the latter case, / again varies linearly with 
slope k until x achieves an extremum. When the mass reverses direction, the contact point again 
sticks and the process repeats in the opposite direction.. 

Given this behavior, the nonlinear force can therefore be expressed as a function of the 
mass's displacement during a cycle of oscillation. The Fourier coefficients can then be calculated 
in terms of the amplitude of vibration. To do this, the displacement is assumed to be sinusoidal 
plus a D.C. offset, i.e. x = Acos 9 + B. Based on this, five different cases can be considered. 
Results for these cases are summarized in the following sections. 

Case 1. No slip : 0 < A < — , and A - — < B < — -A 

k , k k 

d d d 

For this case, an unknown slip may occur during the transient response which results in a 
permanent offset at steady-state. This causes an equally unknown shift in the friction force 
indicated by 7 shown in Figure 2. The D.C. offset B and the mean value of / given by the 
Fourier coefficient f are related by f b = kjl - 7. While they are not known explicitly ( their 
values depend on unknown initial conditions ), they are limited in range, i.e. 

k d A - M N Z jiN- k^A (27) 

In addition, the nonlinear element acts as a spring, i.e. if f » / + / cos 9 + f sin 9 then 

1 - k.A 

c d 

f = o 

s 

and the nonlinear force is in phase with the displacement. 

Case 2. Pure slip : A+B < & +— , A~B £ <5 +— , and A £ — 

1 k 2 k k 

d d d 

This is the case when both the stop distances S Jf & 2 are sufficiently large that slip is not 


47 


hindered. From Figure 3, the friction force 
displacement cycle. 

r ( pN - k^A) + k^Acos 6 

- m N 

f n = < 

(k^A- pN) + k d Acos0 
, h N 


can be divided into four sections during a 


, 0 < 6 < 6 ] 

, 6^ ^ 6 < IT 

(28) 

, n < 0 < 6 2 

, 0 , < e < in 


Where 


0 ( = cos 'Ml' 

d 2 Sn + 6 l 


2 -lH } 


kA 

a 


(29) 


It is observed that the D.C. offset B has no effect on the nonlinear force f in this case. 

n 

Using equations 27 and 28, the nonlinear force f can be expanded in a Fourier series as 

n 

1 * f. + f cos 0 + f sin 0 (30) 

n b c s 


Where 


' b s ° 


/ = - { ( 2 m N - k d A ) sin 0 } + Ar/0, } 

uN uN 
f » - 4 (— — )( 1 - — ) 

5 M 

d 

Case 3. Slip with upper Limit : A + B t 5 + — , £ A- B £ 3 + — 

1 k k 1 2 k 

d d d 


(31) 


In this case (Figure 4), the upper stop is hit during the displacement cycle. Again, the 


nonlinear force can be separated in five sections, as shown below 


f * 

n 


ki B - ) + k Ac os 6 

d id 

- fiN 

( k a A - fiN ) + k^Acos 6 

m N 

i k\B- d.) + k Acos d 

d i d 


, o s e < 0, 

, 0, Z 0 < n 
, n <, 0 < 0 2 
, 0 2 <> 0 K 0 i 
, 0 3 <S 0 < 2n 


( 32 ) 


Where 
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kjl r B)-pN 
8 , = cos { } 


kA 
-] , 2fiN 


0, - 2n ~ cos ' { 1 } 

k.A 

d 

*,( b x ~ B ) + fiN 


8 3 = lit - cos *' { 


k A 

d 


} 


The Fourier coefficients can be derived from equations 32 , 33 as 
1 

f b = —i o ] X 2n + 8 - 9) + pN(9. + 9-26, ) 
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+ k a A( sin 0, + sin 8 , - sin 8 % + 8, - v )} 

1 

f = — f * <5- 3.K sin 0, ~ sin 0, > + k A( 77 + 8. + 8 - 8 + sin 8 , ) 

c 2^ a i j 3d 123 2 

+ fiN( sin 8j + sin 8 3 - 2sin d 2 )} 

uN uN 

f=~2 (—) { + 1 - — } 

* A kA 

d 


Case 4. Slip with Lower Limit : — -S, £ A + B £ 5 + A - B £ 5 + — 

k d ‘ ' k d 2 k d 

In this case, the lower stop is hit during the displacement cycle (Figure 5). As 
nonlinear force f can be separated into four sections. 


f (/i/V- k d A) + * d A:os 8 
- M N 


*,( B + b 2 ) + Ar d -4cos 8 


Where 


\ fiN 


, 2///V 

0 * cos {1 } 

M 

0 


0 2 = 2ff - cos “ 


fiN + Ar (5 + 5 ) 


kA 

d 


} 


8 % = 2ti - cos * ] { 


?N-k(B+ ») 


kA 

d 


1 


0 £ 8 < 8 t 
8 } Z 8 < 8 2 


8 2 Z 8 < 8 3 

8 3 £ 8 £ 2tt 


( 33 ) 


(34) 


before, the 


(35) 


( 36 ) 


From equation 35 and 36, the Fourier coefficients are obtained. 
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f . = { pNl 28 ~ 0, - 6 + 2v ) + k A( sin 8 + sin 8 ~ sin 0 - 8 ) 


2n 


k ( B + a,)< 8-8,)} 


f'. = — { pN( 2sin 8 ] - sin 8 , - sin 8 % ) + *,( B + <5, )( sin 0 3 - sin 8 , ) (37) 

2 77 

+ k iJ A(8 ] + 8 3 -8 2 -sin8 t ) } 

uN B+ *2 uN 

f = -2 <— >{ — + 1 } 

s v A k A 

d 


Case 5. Slip with Both Limit : A + B £ , A - B £ + — 


In this case, both stops are hit during the displacement cycle( Figure 6). As before, the 


nonlinear force i can be divided into five sections. 

n 


f 

n 


Where 


f k ( B - d. ) + k.Acos 8 

d 1 d 

- M N 

< k ( B + d , ) + k >4cos 8 

d 1 d 

ftN 

i kj B- 5 j ) + * d 4cos 8 


8 cos { 


k d (d-B)- M H 


kA 

d 


8 2 - cos ■' { - 


8 , - 2n - cos { 


8 = 2tr ~ cos ** { 

4 


pN + k ( B + i ) 


kA 

d 


} 


ftN- k <B + 6 ) 


kA 

d 


-> 


HN + k ( & - B) 


kA 

d 


-} 


, 0 < 8 < 0 , 

, 8 i Z 8 < 8 2 

, 8 2 < 8 < e 3 
, 8 3 £ 8 < 8 a 
, 0 4 Z 8 < 2n 


(38) 


( 39 ) 


And the Fourier coefficients are derived from 38, 39 as 
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where 


f= — { ki b- s, >< e , + 2n - e j + k \ b+ o)( e - 0 , > 

b«d tl 4 d 232 

2v 

+ k^Ai sin 8 t + sin 6 3 ~ sin 0, - sin 0 4 ) 

+ fiN{ e t + 0 4 - e 2 - e 3 » ' 


f - — { k < B - 3. )( sin 8 - sin 0, ) + k A( 8+ 8+2n - 8 - 6 ) 

c ~ <J J 1 4 d 1 3 24 

2n 

+ k ( B + 3 0 )( sin 0 - sin 0 ) 

0 2 3 2 

+ fiNi sin d ] - sin 0 ^ + sin 0 4 - sin 0 3 ) } 

2 <5 

f = (— ) uN where 5 = 3 + 

s * ^ 


From the above obtained results, the nonlinear force f can be written as 

n 

f Ti =f h (A,B' i iN'k^&^h 2 ) + f{A'B y tiN > k^b^b 2 )cos 0 
+ f{ A , B , /i/V , /r d , , 3^ )sin 0 

The three coefficients can be nondimensionhzed as follows : 

f b m fNF b (a,fi,y.q) 

f c - F c < A r> <7> 

f = A ( a ' A y* 

k A A 

a 

a ~ 

fiN 


B 



where, 3 = ^ 1 + ^ 2 


MO) 


(41) 


(42) 


( 43 ) 


